home *** CD-ROM | disk | FTP | other *** search
/ Video Toaster 4.0 / Video Toaster v4.0.iso / arexx / modeler / plot2d.lwm < prev    next >
Text File  |  1993-12-13  |  12KB  |  474 lines
  1. /* CMD: Surface Plot
  2.  *
  3.  * 2-D Function Surface Maker for Modeler
  4.  * Originally by Arnie Cachelin © 1992 NewTek Inc., Sat Jul 11 1992
  5.  *
  6.  * Modified to use Modeler-based UI by Stuart Ferguson, 11/92
  7.  * Further modified by Arnie, Sat May  8 16:04:25 1993
  8.  */
  9.  
  10. arg bas
  11. call addlib "LWModelerARexx.port", 0
  12. call addlib "rexxsupport.library", 0, -30, 0
  13. signal on error
  14. signal on syntax
  15.  
  16. MATHLIB="rexxmathlib.library"
  17. IF POS(MATHLIB , SHOW('L')) = 0 THEN
  18.   IF ~ADDLIB(MATHLIB , 0 , -30 , 0) THEN DO
  19.     call notify(1,"!Can't find "MATHLIB)
  20.     exit
  21.     END
  22.  
  23. sysnam = 'Plot 2D Function'
  24. filnam = 'ENV:plot2d.state'
  25. version = 'Plot 2D v1.1'
  26.  
  27. /* Setup state.  Read stored one, if any.
  28.  */
  29. NSX  = 20
  30. NSY  = 20
  31. xmin = -10
  32. xmax = 10
  33. ymin = -10
  34. ymax = 10
  35. func = '3*sin(x*y)'
  36. ufunc=func
  37. fn = 3
  38. flip = 1
  39. tri = 1
  40. if bas="" then bas = 0
  41. bas = 0
  42. typ=2
  43. if (exists(filnam)) then do
  44.     if (~open(state, filnam, 'R')) then break
  45.     if (readln(state) ~= version) then break
  46.     parse value readln(state) with nsx nsy xmin xmax ymin ymax fn tri typ .
  47.     func = readln(state)
  48.     call close state
  49. end
  50.  
  51. FnList.1='RandomSheet'
  52. FnList.2='Wave'
  53. FnList.3='RadialWave'
  54. FnList.4='Gaussian'
  55. FnList.5='Interfere'
  56. FnList.6='Custom'
  57. FList= FnList.1 FnList.2 FnList.3 FnList.4 FnList.5 FnList.6
  58.  
  59. call req_begin sysnam
  60.  
  61. id_fnl = req_addcontrol("f(x,y)", 'CH',FList)
  62. id_lox = req_addcontrol("Low X", 'n', 1)
  63. id_hix = req_addcontrol("High X", 'n', 1)
  64. id_nsx = req_addcontrol("X Segments", 'n')
  65. id_loy = req_addcontrol("Low Y", 'n', 1)
  66. id_hiy = req_addcontrol("High Y", 'n', 1)
  67. id_nsy = req_addcontrol("Y Segments", 'n')
  68. id_fun = req_addcontrol("Custom Function", 's', 35)
  69. id_tri = req_addcontrol("Triangles", 'b')
  70. id_typ = req_addcontrol("Build: ","CH","Points Polys Curves")
  71.  
  72. call req_setval id_lox, xmin, -10
  73. call req_setval id_loy, ymin, -10
  74. call req_setval id_hix, xmax, 10
  75. call req_setval id_hiy, ymax, 10
  76. call req_setval id_nsx, nsx, 20
  77. call req_setval id_nsy, nsy, 20
  78. call req_setval id_fun, func
  79. call req_setval id_fnl, fn,6
  80. call req_setval id_tri, tri,tri
  81. call req_setval id_typ, typ,typ
  82.  
  83. if (~req_post()) then do
  84.     call req_end
  85.     exit
  86. end
  87.  
  88. NSX = req_getval(id_nsx) % 1
  89. NSY = req_getval(id_nsx) % 1
  90. xmin = req_getval(id_lox)
  91. xmax = req_getval(id_hix)
  92. ymin = req_getval(id_loy)
  93. ymax = req_getval(id_hiy)
  94. fn   = req_getval(id_fnl)
  95. typ   = req_getval(id_typ)
  96. if fn=6 then func = req_getval(id_fun)
  97. else func = FnList.fn'(x,y)'
  98. tri = req_getval(id_tri)
  99.  
  100. call req_end
  101.  
  102. if (open(state, filnam, 'W')) then do
  103.     call writeln state, version
  104.     call writeln state, nsx nsy xmin xmax ymin ymax fn tri typ
  105.     call writeln state, ufunc
  106.     call close state
  107. end
  108.  
  109. xrange = xmax - xmin
  110. yrange = ymax - ymin
  111. xmesh = xrange / NSX
  112. ymesh = yrange / NSY
  113. rscale = sqrt(xrange*yrange)
  114. tri_height = sqrt(3)/2
  115. ifunc = "z =" func
  116. say ifunc
  117. z=0
  118. zmax=z
  119. zmin=z
  120. call randu(time('s'))  /* Seed random number generator */
  121.  
  122. if typ=3 then call RectCurves
  123. else
  124. if tri then do
  125.   if typ=1 then call TriPoints
  126.   else call TriMesh
  127.   end
  128. else do
  129.   if typ=1 then call RectPoints
  130.   else call RectMesh
  131.   end
  132.  
  133. /* if Bas totalpoints=totalpoints+MakeBase() */
  134.  
  135. l1 = "Points created:" totalpoints
  136. l2 = "Polygons created:" poly
  137. l3 = "Z ranges between" zmin "and" zmax
  138. call notify 1, '!'sysnam, l1, l2, l3
  139.  
  140. exit
  141.  
  142. RectMesh:
  143.     totalpoints = (NSX+1) * (NSY+1)
  144.     totalpolys  = NSX * NSY
  145.     call add_begin
  146.     call meter_begin totalpoints+2, sysnam, "Computing "totalpoints" points for "totalpolys" squares"
  147.     do y=ymin to ymax by ymesh
  148.       if y=ymax then TopCorner.4=totalpoints
  149.       do x=xmin to xmax by xmesh
  150.         interpret ifunc
  151.         if z<zmin then zmin=z
  152.         if z>zmax then zmax=z     /* Just some silly stats for later */
  153.         if (flip) then vec =x z y
  154.         else vec = x y z
  155.         call add_point(vec)
  156.         call meter_step
  157.       end
  158.       if y=ymin then TopCorner.2=totalpoints
  159.       if y=ymax then TopCorner.3=totalpoints
  160.     end
  161.  
  162.     point=1
  163.     poly=0
  164.     call meter_begin totalpolys, sysnam, "Generating "totalpolys" Polygon Mesh"
  165.     do y=ymin to ymax-ymesh by ymesh   /* Don't make wrap-around polygon */
  166.       do x=xmin to xmax by xmesh
  167.         if x<xmax then do   /* Again, Don't make wrap-around polygon! */
  168.           if (flip) then
  169.             call add_quad point point+NSX+1 point+NSX+2 point+1
  170.           else
  171.             call add_quad point point+1 point+NSX+2 point+NSX+1
  172.           poly = poly + 1
  173.           call meter_step
  174.         end
  175.         point = point + 1
  176.       end
  177.     end
  178.     call meter_end
  179.     call add_end
  180. return totalpoints
  181. /*  */
  182.  
  183. RectPoints:
  184.     totalpoints = (NSX+1) * (NSY+1)
  185.     poly=totalpoints
  186.     call add_begin
  187.     call meter_begin totalpoints+2, sysnam, "Computing "totalpoints" points for 1-point polygons"
  188.     point=1
  189.     do y=ymin to ymax by ymesh
  190.       do x=xmin to xmax by xmesh
  191.         interpret ifunc
  192.         if z<zmin then zmin=z
  193.         if z>zmax then zmax=z     /* Just some silly stats for later */
  194.         if (flip) then vec =x z y
  195.         else vec = x y z
  196.         call add_point(vec)
  197.         call add_polygon(point)
  198.         point=point+1
  199.         call meter_step
  200.       end
  201.     end
  202.     call add_end
  203. return totalpoints
  204. /*  */
  205.  
  206. RectCurves:
  207.     totalpoints = (NSX+1) * (NSY+1)
  208.     poly=NSX+1 + NSY+1
  209.     call add_begin
  210.     call meter_begin totalpoints+NSX+2, sysnam, "Computing "totalpoints" points for "poly" Curves"
  211.     crv=""
  212.     point=1
  213.     do y=ymin to ymax by ymesh
  214.       do x=xmin to xmax by xmesh
  215.  
  216.         interpret ifunc
  217.  
  218.         if z<zmin then zmin=z
  219.         if z>zmax then zmax=z
  220.  
  221.         if (flip) then vec =x z y
  222.         else vec = x y z
  223.         call add_point(vec)
  224.         crv=crv point
  225.         point=point+1
  226.         call meter_step
  227.       end
  228.       call Add_Curve(crv)
  229.       crv=""
  230.     end
  231.     do p=1 to NSX+1
  232.       do o=0 to NSY
  233.         crv=crv p+o*(NSX+1)
  234.       end
  235.       call meter_step
  236.       call Add_Curve(crv)
  237.       crv=""
  238.     end
  239.     call add_end
  240. return totalpoints
  241. /*  */
  242.  
  243. TriMesh:
  244.     totalpoints = (NSX+1) * (NSY+1)
  245.     totalpolys  = NSX * NSY * 2
  246.     call add_begin
  247.     call meter_begin totalpoints*2, sysnam, "Computing "totalpoints" points"
  248.     offset=0
  249.     rows=0
  250.     totalpoints=0
  251.  
  252.     do y=ymin to ymax by ymesh* tri_height
  253.         rows=rows+1
  254.         columns=0
  255.       if y=ymax then TopCorner.4=totalpoints
  256.       do x=xmin+offset to xmax+offset by xmesh
  257.  
  258.              columns = columns + 1
  259.         interpret ifunc
  260.  
  261.         if z<zmin then zmin=z
  262.         if z>zmax then zmax=z     /* Just some silly stats for later */
  263.  
  264.         if (flip) then vec =x z y
  265.         else vec = x y z
  266.         call add_point vec
  267.         totalpoints = totalpoints + 1
  268.         call meter_step
  269.       end
  270.       if y=ymin then TopCorner.2=totalpoints
  271.       if y=ymax then TopCorner.3=totalpoints
  272.         if offset=0 then offset=.5 * xmesh  /* offset alternate lines */
  273.         else offset=0
  274.     end
  275.     call meter_end
  276.  
  277.     point=1
  278.     poly=0
  279.     off=0
  280.     call meter_begin totalpolys, sysnam, "Generating "totalpolys" Polygon Mesh"
  281.     do row=0 to rows-2
  282.       if off=0 then off=1  /* Boy this feels kludgey!!! */
  283.       else off=0
  284.          do col=1 to columns - 1
  285.         if (flip) then do
  286.           call add_quad col+row*columns col+(row*columns)+1 col+((row+1)*columns)+abs(off-1)
  287.           call add_quad col+(row*columns)+off col+((row+1)*columns)+1 col+((row+1)*columns)
  288.           poly=poly+2
  289.           end
  290.         else do
  291.           call add_quad col+row*columns col+((row+1)*columns)+abs(off-1) col+(row*columns)+1
  292.           call add_quad col+(row*columns)+off col+((row+1)*columns) col+((row+1)*columns)+1
  293.           poly=poly+2
  294.           end
  295.         call meter_step
  296.       end
  297.     end
  298.     call meter_end
  299.     call add_end
  300. return totalpoints
  301. /*  */
  302.  
  303. TriPoints:
  304.     totalpoints = (NSX+1) * (NSY+1)
  305.     call add_begin
  306.     call meter_begin totalpoints*2, sysnam, "Computing "totalpoints" points for 1-point polygons"
  307.     offset=0
  308.     rows=0
  309.     totalpoints=0
  310.     point=1
  311.  
  312.     do y=ymin to ymax by ymesh* tri_height
  313.         rows=rows+1
  314.         columns=0
  315.       do x=xmin+offset to xmax+offset by xmesh
  316.  
  317.              columns = columns + 1
  318.         interpret ifunc
  319.  
  320.         if z<zmin then zmin=z
  321.         if z>zmax then zmax=z     /* Just some silly stats for later */
  322.  
  323.         if (flip) then vec =x z y
  324.         else vec = x y z
  325.         call add_point vec
  326.         call add_polygon point
  327.         point=point+1
  328.         call meter_step
  329.       end
  330.         if offset=0 then offset=.5 * xmesh  /* offset alternate lines */
  331.         else offset=0
  332.     end
  333.     call meter_end
  334.     totalpoints = point - 1
  335.     poly=totalpoints
  336.     call add_end()
  337. return totalpoints
  338. /*  */
  339.  
  340.  
  341. MakeBase:
  342.     add_begin
  343.     point=TotalPoints-1
  344.     say point
  345.     TopCorner.1=1
  346.     BotCorner.1=totalpoints
  347.     BotCorner.2=BotCorner.1 + NSX +1
  348.     BotCorner.3=BotCorner.2 + NSY
  349.     BotCorner.4=BotCorner.3 + NSX +1
  350.     z=ZMin-(xrange+yrange)/4
  351.     y=ymin
  352.     BotCorner.1=point+1
  353.     do x=xmin to xmax by xmesh
  354.       call add_point x y z
  355.          point=point+1
  356.        pointList.point = x y z
  357.     end
  358.     BotCorner.2=point
  359.     x=xmin
  360.     do y=ymin to ymax by ymesh
  361.       call add_point x y z
  362.          point=point+1
  363.        pointList.point = x y z
  364.     end
  365.     y=ymax
  366.     BotCorner.4=point+1
  367.     do x=xmin to xmax by xmesh
  368.       call add_point x y z
  369.          point=point+1
  370.        pointList.point = x y z
  371.     end
  372.     BotCorner.3=point
  373.     x=xmax
  374.     do y=ymin to ymax by ymesh
  375.       call add_point x y z
  376.          point=point+1
  377.        pointList.point = x y z
  378.     end
  379.  
  380.     call surface("AreaBottom")
  381.     call add_quad BotCorner.1 BotCorner.2 BotCorner.3 BotCorner.4
  382.     poly=poly+1
  383.     call surface("AreaBase")
  384.  
  385.     do i=0 to NSX-1
  386.       call add_quad TopCorner.1+i TopCorner.1+i+1 BotCorner.1+i+1 BotCorner.1+i
  387.       poly=poly+1
  388.     end
  389.  
  390.     do i=0 to NSY-1
  391.       call add_quad TopCorner.4+i TopCorner.4+i+1 BotCorner.4+i+1 BotCorner.4+i
  392.       poly=poly+1
  393.     end
  394.  
  395.     do i=0 to NSX-1
  396.       a=((TopCorner.1)+i*(NSX+1))
  397.         b=((TopCorner.1)+(i+1)*(NSX+1))
  398.         c=((BotCorner.2)+i+1)
  399.         d=((BotCorner.2)+i+2)
  400.       call add_quad a b d c
  401.       poly=poly+1
  402.     end
  403.  
  404.     do i=0 to NSY-1
  405.       a=((TopCorner.2)+i*(NSY+1))
  406.         b=((TopCorner.2)+(i+1)*(NSY+1))
  407.         c=((BotCorner.3)+i+1)
  408.         d=((BotCorner.3)+i+2)
  409.       call add_quad a b d c
  410.       poly=poly+1
  411.     end
  412.  
  413.     add_end
  414. return point
  415.     /*  */
  416.  
  417. Radius: PROCEDURE
  418.   arg xf, yf
  419.   return sqrt(xf*xf+yf*yf)
  420.  
  421. Sinc:   PROCEDURE EXPOSE rscale
  422.   arg xf, yf     /* Classic the spherical Bessel f'n j0 */
  423.   r=Radius(xf,yf)
  424.   if (r = 0)
  425.     zf = rscale
  426.   else
  427.     zf = rscale * sin(r) / r
  428.   return zf
  429.  
  430. Wave:   PROCEDURE  EXPOSE rscale
  431.   arg xf, yf     /* Simple wavy sheet */
  432.   zf=rscale*sin(2*3.141592*xf/rscale)/3
  433.   return(zf)
  434.  
  435. RadialWave:   PROCEDURE  EXPOSE rscale
  436.   arg xf, yf     /* rings of waves sheet */
  437.   xc=0; yc=0     /* Center coord.s */
  438.   r=Radius(xf-xc,yf-yc)
  439.   zf=rscale*sin(10*3.141592*r/rscale)/5
  440.   return(zf)
  441.  
  442. Interfere:   PROCEDURE  EXPOSE rscale
  443.   arg xf, yf     /* Interference of several (3) sources */
  444.   r=Radius(xf,yf)      /* (0,0)  */
  445.   d1=Radius(xf-5*rscale/10,yf+2*rscale/10) /* (5,-2) */
  446.   d2=Radius(xf+4*rscale/10,yf-3*rscale/10) /* (-4,3) */
  447.   zf=(rscale/7)*(sin(r)+sin(d1*1.5)+1.4*cos(d2))
  448.   return(zf)
  449.  
  450. RandomSheet:   PROCEDURE  EXPOSE rscale
  451.   arg xf, yf     /* Random altitudes */
  452.   amp=(rscale/3)
  453.   zf=randu()*amp-amp/2
  454.   return(zf)
  455.  
  456. Polynomial:   PROCEDURE  EXPOSE rscale
  457.   arg xf, yf
  458.   zf = 2*xf*xf + 3*xf - 2*yf*yf + 8*yf + 2
  459.   zf = zf / 20
  460.   return(zf)
  461.  
  462. Gaussian:   PROCEDURE  EXPOSE rscale
  463.   arg xf, yf
  464.   xc=0; yc=0     /* Center coord.s */
  465.   amp = rscale/3
  466.   r = Radius(xf-xc,yf-yc)
  467.   zf = amp * exp(-r*r/20)
  468.   return(zf)
  469.  
  470. syntax:
  471. error:
  472.   call end_all
  473.     t=Notify(1,'!Rexx Script Error','@'ErrorText(rc),'Line 'SIGL)
  474.     exit